A dynamical mean-field theory study of Nagaoka ferromagnetism. 
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We revisit Nagaoka ferromagnetism in the [/ = oo Hubbard model within the dynamical mean- 
field theory (DMFT) using the recently developed continuous time quantum Monte Carlo method 
as the impurity solver. The stability of Nagaoka ferromagnetism is studied as a function of the 
temperature, the doping level, and the next-nearest-neighbor lattice hopping t' . We found that the 
nature of the phase transition as well as the stability of the ferromagnetic state is very sensitive 
to the t' hopping. Negative t' = —Q.lt stabilizes ferromagnetism up to higher doping levels. The 
paramagnetic state is reached through a first order phase transition. Alternatively, a second order 
phase transition is observed at t' — 0. Very near half-filling, the coherence temperature T^oh of the 
paramagnetic metal becomes very low and ferromagnetism evolves out of an incoherent metal rather 
than conventional Fermi liquid. We use the DMFT results to benchmark slave-boson method which 
might be useful in more complicated geometries. 



I. INTRODUCTION 

The stability of the ferromagnetic phase in the U = oo 
Hubbard model is a long standing problem. Nagaokai 
showed that for a single hole in a bipartite lattice the 
ground state is a fully polarized fcrromagnct, and the 
term " Nagaoka ferromagnetism" is commonly used to de- 
scribe this state. Whether a fully or a partially polarized 
phase persist to a finite hole density (S) is controversial 
and has been the subject of numerous investigations^. 

The problem has been addressed with variational wave 



functions 
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slave particle methods — i^, quantum 



Monte Carlo (QMC) methodsiS, and variational QMC 
methodaii. In all these methods the ferromagnetism 
is stable up to a critical value of doping 5c- It was 
also demonstrated by these approaches that the size of 
the ferromagnetic region depends strongly on the lat- 
tice through the electronic dispersion. The ferromagnetic 
state was found to be unstable even for the case of a sin- 
gle hole in the U = oo square lattice with a small positive 
next-nearest neighbor hopping t'—. At an intermediate 
or a large U, a flat band below the Fermi level— or a peak 
in the density of states below the Fermi level I4,i5^i6,i7 ^ 
realized in the fee lattic o^^i^^ or a Van Hove singularitj«2£, 
stabilize the ferromagnetic state. 

The dynamical mean-field theory (DMFT) has also 
been used to address the Nagaoka problem, however the 
number of available impurity solvers in the U = oo case 
is very limited. Obermeier et al^ carried out the first 
DMFT study of this problem using the non-crossing ap- 
proximation as the impurity solver. They found a par- 
tially polarized ferromagnetic state below a critical tem- 
perature Tc in the infinite dimensional hypercubic lattice. 
The existence of a ferromagnetic state in this model was 
later confirmed by a DMFT study which used numerical 
renormalization group as the impurity solver—. 

In this study, we revisit the problem of Nagaoka fer- 
romagnetism in the U = oo Hubbard model within 
DMFT, using the recently developed continuous time 
quantum Monte Carlo (CTQMC) method as the impu- 



rity solveri22i2i. This impurity solver allows the numeri- 
cally exact solution of the DMFT equations at very low 
temperatures for all values of doping level 5 even in the 
U = oo model. We find that at large doping, the fer- 
romagnetism emerges from a conventional Fermi liquid, 
while at small doping the Curie temperature is very close 
to the coherence temperature, hence the ferromagnetism 
emerges from an incoherent state. We pay particular 
attention to the possibility of phase separation and its 
dependence on the sign of t' /t. Finally we benchmark 
simpler approaches to the problem such as the slave bo- 
son method. Within slave boson approach, several phys- 
ical quantities such as the quasiparticle renormalization 
amplitude or the susceptibility can not be determined 
reliably. Nevertheless we show that the total energy can 
be computed quite reliably within the simple slave bo- 
son approach due to error cancellation. This is impor- 
tant since the detailed modeling of optical lattices of cold 
atoms, which provide a clean realization of the Hubbard 
model, will require incorporating spatial inhomogeneities 
into the treatments of strong correlations. At present, 
this can only be done with simpler techniques such as 
slave bosons methods. 

We study the Hamiltonian of the U ^ oo Hubbard 
model given by 



(1) 



where Pg is a projection operator which removes states 
with double-occupied sites. We choose the lattice dis- 
persion of the two dimensional square lattice with the 
nearest-neighbor (n.n) hopping t and the next-nearest- 
neighbor (n.n.n) hopping t'. The units are fixed by choos- 
ing t = i. 



II. A DMFT-fCTQMC APPROACH 

DMFT maps the partition function of the Hubbard 
model onto the partition function of an effective Ander- 
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son impurity model (AIM) resulting in the following ef- 
fective action. 

Seff = Satom+ cLt dr' ^ 4 (t ) A, (t - t')c, (t' )(2) 

Jo Jo ^ 

where Satom represents the action of the isolated impu- 
rity, and Act(t — t') is the hybridization hmction of the 
effective AIM. In this U ^ oo case, the double occupied 
state of the impurity should be excluded when evaluat- 
ing Satom- Ao-(t — t') is uot initially known and it must 
be determined by the DMFT self-consistency condition 
given below. The impurity Green function and the im- 
purity self-energy are given by the following equations 

GAr-r') = -{Tc,ir)cUr'))s.„ (3) 

T,aiiuJn) = iljJn + M ^ Acr(icj„) - ^^^{iuJn)- (4) 

The DMFT self-consistency condition requires that the 
local Green's function of the lattice coincides with the 
Green's function of the auxiliary AIM and identifies the 
equivalence between the lattice local self-energy and the 
self-energy of the corresponding AIM, i.e., 

y I 

^ ILUn + /-t + /icr - e(k) - Y.a{lUJn) 



iuJn + f-l + ha - Aa{iuJn) - Ea(icJ„) ' 

where the lattice dispersion of our choice is e(k) = 
—2t{coskx + cosfcy) — 4t' cos fca; cos fcy and h is the ex- 
ternal magnetic field. For a given Weiss field Acr(iu;„), 
the effective action Se/f is constructed and the AIM is 
solved for the new Ga{iuJn) and S]cr(«^ra)- Using the self- 
consistency condition Eq.5, the new Weiss field /\„{iojn) 
is computed. This iterative procedure is repeated until 
the Green's function is converged. 

To solve the impurity problem of Eq. d the CTQMC 
impurity solver is used. In this method, the hybridization 
part of the effective action is treated as a perturbation 
around the atomic action and all diagrams are summed 
up by stochastic Metropolis samplingi^ In this U = oo 
case, doubly occupied state of the atom is excluded from 
atomic eigenstates. CTQMC converges well in the low 
Matsubara frequency region, but it is poorly behaved in 
the high frequency region. Therefore, one needs the ana- 
lytic expression for the self-energy in the high frequency 
limit and it has to be interpolated to the low frequency 
region. The high frequency expansion for the U ~ oo 
Hubbard model gives 

Re[Y.^{oo)] = mi^/ml^ + /i (6) 
Im[Y.aioo)] = {I - l/mo^)uj (7) 

where moa = ({C(t,4}) = 1 - "--(t, rrii^ = 
{{[c„, H],cl}) = - n_,) - Tr[A_,G_,]. Note the 

appearance of the kinetic energy Tr[A^^G-a] in this ex- 
pansion which is absent for finite U. 



Within CTQMC; various spin dependent physical 
quantities can be calculated such as occupation numbers 
(n|,n|) and the local magnetic susceptibility (xioc)- The 
q ~ magnetic susceptibility of a lattice can be cal- 
culated from xioc by evaluating the two particle vertex 
functions, which is a numerically demanding task. To cir- 
cumvent this difficulty, Xq=o of a lattice can be calculated 
from the ratio of magnetization to the external magnetic 
field (x = ^\h=o)- The external field h alters the ef- 
fective action (Eq. [J) by adding ha to atomic energies 
and the self-consistency condition (Eq. [5]) is enforced to 
include the spin dependent ha term during DMFT iter- 
ations. The exclusion of the double occupancy (C/ = oo) 
implies the Hubbard potential energy to vanish and the 
only relevant energy is the kinetic energy. The latter is 
given by Tr[Aa-Gcr], and it is related to the average of 
the perturbation order k as follows: 

EMn,a - Tr[A,G^] = -T{k,) (8) 

where T is temperature. Therefore, it is possible to cal- 
culate the kinetic energy to high accuracy by evaluating 
{k„). The free energy, _F, can also be derived from the 
kinetic energy as long as the system is in the Fermi liquid 
regime 

F{T) ^ Ek,n - y Z-Vo(A')r', (9) 

where Z is the renormalization residue and po is the non- 
interacting density of states. 

Fig.[TJ(a) shows the reduced magnetization 'mr~{n-^ — 
ri|)/(ri,| + ni) as a function of the electron density n at 
three distinct t' /t ratios. The result is notably differ- 
ent for different values of t' /t. The spontaneously bro- 
ken ferromagnetic (FM) state {rrir ^ 0) is favored for 
t'/t < while the FM state is unstable for t' /t > 0. The 
critical density (ric) at which the transition occurs in- 
creases as t'/t increases, reducing the region of stability 
of the FM state. Moreover, at t'/t = —0.1 magnetization 
rrir changes abruptly at ric=0.705 indicating a first order 
transition, while at t'/t = magnetization rrir increases 
continuously indicating a second order phase transition 
at nc=0.815. 

Notice that close to half filling the Curie temperature 
is low and at fixed temperature (T = 0.01) it becomes in- 
creasingly difficult to converge the DMFT equations near 
the transition temperature due to the standard critical 
slowing down. 

Near half filling the quasiparticle bandwidth is small 
due to strong correlations hence the thermal fluctuations 
are comparable to the Curie temperature in this region. 
A stable FM state is possible only if T is sufficiently 
lower than Tcoh- In the region above 0.95, an incoherent 
paramagnetic (PM) state becomes stable as T exceeds 

J^coh- 

Inspecting the chemical potential as a function of den- 
sity reveals that the nature of the transition changes 
with t'/t (see figure [T](b)). For t'/t = the transition 
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FIG. 1: (Color online) (a) The reduced magnetization 
mr={ni — n^)/{ni + n^) vs the electron density n at t' /t=- 
0.1, 0, and 0.1 (b) the chemical potential fj, vs n at t' /t=-0.1, 
0, and 0.1. Filled points indicate a FM state. Inset : FM 
free energy and PM free energy vs n at t' /t—-Q.l. The dot- 
ted line is constructed using the Maxwell construction. All 
calculations were performed at T=0.01. 



is continuous while for t'/t = —0.1, there is a region 
of constant chemical potential which corresponds to a 
first-order transition. The flat chemical potential region 
(n = 0.696 - 0.715) indicates that two different DMFT 
solutions (FM. PM) can be converged depending on the 
initial conditions and it indicates phase separation (PS) 
of the FM and PM state. This region is determined by 
Maxwell construction which connects common tangents 
between two phases in the free energy vs n graph. (Fig. [1] 
inset) 

The original debate on the Nagaoka problem was fo- 
cused on the existence of the fully polarized FM state at 
finite S in the T — )■ limit. Therefore, it is necessary to 
investigate magnetization at very low T. In Fig. [2] we 
show very low temperatures (T = O.OOK) results and it is 
clear that the magnetization saturates to a value smaller 
than unity for t'/t = while it reaches unity at low tem- 
peratures for t'/t = —0.1. The fully polarized Nagaoka 
state is thus not stable for t'/t — and moderately small 
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FIG. 2: (Color online) m,. vs T at fixed n =0.85 with t' /t^- 
0.1 and 0. The fully polarized FM state (m^ — 1) is expected 
only when t' /t=-0.1. 



doping {S ^ 0.1) while it is realized for t'/t = —0.1. As 
the spins become fully polarized (t'/t = — 0.1,r 0), 
numerics requires high statistics and an error-bar is spec- 
ified to take into account the numerical error. 

The spectral functions are shown in Fig. [31 Since 
CTQMC delivers response functions on the imaginary 
frequency axis, one needs to perform the analytical con- 
tinuation of the Green function to the real axis. Here 
we use the maximum entropy method^^. The spectral 
functions show noticeable differences for small change in 
t'. At t'/t = —0.1, the majority spin spectral function 
shows a very small renormalization due to interactions 
(Z ~ 1) and a large spectral peak in the occupied part 
of the spectra. The overall shape is similar to the non- 
interacting spectral function (Fig. [31 inset). The minor- 
ity spin spectral function is much more correlated and 
shows a narrow quasiparticle band above the Fermi level 
and a tiny lower Hubbard band. In the magnetic state, 
the occupied part of the spectra is thus well described by 
a model of a weakly correlated FM metal. 

At t'/t ~ and t'/t = 0.1, the spectral functions con- 
sist of both the narrow quasiparticle band and the lower 
Hubbard band. In the U = oo Hubbard model, the upper 
Hubbard band disappears due to the exclusion of double 
occupancy. 

The stabihty of the FM state at t'/t = -0.1 can be 
traced back to the large spectral peak in the occupied 
part of the spectra of the non-interacting DOS shown in 
the inset of Fig. [31 As explained above, the majority spin 
of the FM state shows only weak renormalization due to 
interactions. This is a consequence of the Pauli exclusion 
principle which constrains the motion of a hole in the po- 
larized background and interactions, being less important 
in this case, do not hamper the coherent motion of the 
hole through the polarized background. The kinetic en- 
ergy of this state thus clearly depends on t'/t ratio and is 
reduced with decreasing t'/t. Contrary to the FM state. 
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FIG. 4: (Color online) The uniform susceptibility (Xqlo) 
n at t' /t = and 0.1. The dotted line is for the extrapolation 
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FIG. 5: (Color online) The critical temperature Tc vs n at 
t' = 0. ric at T = is obtained from the extrapolation. The 
dotted line represents the coherence temperature Tcoh vs n. 



FIG. 3: (Color online) The spectral functions A{u}) at t' /t=- 
0.1 (top), (middle), and 0.1 (bottom) for fixed n =0.85. 
Inset: Non- interacting spectral functions {Ao{ui)) of the ma- 
jority spin at the corresponding t' /t values, {fio = /x— _ReE(0)) 
All calculations were performed at T=0.01. 



the correlations are very strong in the PM state regard- 
less of the speetral peak in the non-interacting DOS and 
t' /t ratio. The coherent part of the spectra does not con- 
tribute much to the kinetic energy as the quasiparticle 
bandwidth shrinks due to the strong correlations. The 
incoherent part of the spectra in the form of the Hubbard 
bands arises from localized electrons and consequently it 
is almost independent of the specific lattice dispersion. 
Therefore, the kinetic energy of the PM state weakly de- 



pends on t' /t ratio. The peak in the occupied part of 
the spectra of the non-interacting DOS thus reduces the 
kinetic energy of the FM state compared to the PM state 
thus stabilizing ferromagnetism. 

It is known from other studiesi^ that a highly degen- 
erate flat band in the occupied part of the spectra favors 
ferromagnetism at any finite U. However, this flat band 
ferromagnetism (an extreme limit of the Stoner ferromag- 
netism) argument is not applicable to the t' /t — —0.1 
case of the Nagaoka ferromagnetism (the other extreme 
limit of the Stoner ferromagnetism). In a flat-band 
model, the ground state of the non-interacting system 
is highly degenerate due to the presence of the flat band. 
However, even a small Coulomb repulsion lowers the en- 
ergy of the FM state (if the flat band is half-filled) and 
stabilizes the FM state. The role of the Coulomb in- 
teraction is simply to lift the huge degeneracy and "se- 
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FIG. 6: (Color online) The local susceptibility (Xioc) ^ ^'^'^ 
the uniform susceptibility (Xqlo) vs T (t' jt = 0) 
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lect" the states with the highest magnetization as unique 
ground states. In the infinite \J model, the potential en- 
ergy vanishes because of no doubly occupancy. However, 
the kinetic energy depends sensitively on the smoothness 
of the spin polarized background, and a disordered PM 
state can not gain the kinetic energy by the variation of 
t' jt while a FM state can. 

The inverse of the uniform magnetic susceptibihty 
(x^=o) °f ^'^^ state vs n at t' jt = and 0.1 is shown 
in Fig. m The extrapolated line at t' = indicates that 
X diverges near n = 0.815, confirming the second or- 
der transition at the critical density [n^ = 0.815). At 
t' /t = 0.1, one might expect x will diverge near n = 1. 
However, as Tcoh becomes smaller than T near n = 1, 
the incoherent PM state is stabilized. In other words, at 
t' /t = 0.1, the crossover from the coherent PM state to 
the incoherent PM state occurs instead of the transition 
to the FM state. 

Fig. [5] shows the critical temperature (Tc) vs n at 
t' /t = 0. In the region below Tc a partially polarized 
FM state is found, and it is determined by observing 
n| 7^ in a CTQMC result. This graph shows that 
the lower critical density (ric) at T = is around 0.8. 
At half filling critical temperature should vanish due to 
the following reason: The kinetic energy at half filling is 
zero in both the PM and the FM state because of the 
blocking of charge density. The entropy of the paramag- 
net is much larger than the entropy of the fcrromagnct 



FIG. 7: (Color online) ImT,{iuj„) vs iun at the coherent 
FM state (the top panel), the coherent PM state (the mid- 
dle panel), and the incoherent PM state (the bottom panel). 
{t'/t = 0) 



due to the large spin degeneracy of the PM state. In 
other words, PM state is thermodynamically stable at 
any finite temperature at n = 1. 

As the width of the quasiparticle band becomes smaller 
near ti = 1, the coherence temperature Tcoh is also re- 
duced making it hard to sustain the quasiparticle coher- 
ent band. At T > Tcoh, the PM state is clearly stabilized. 
The Tcoh boundary can be determined from the imagi- 
nary part of self energy {ImT,{iujn)) on the imaginary 
frequency axis. In a coherent region {T < Tcoh), the 
renormalization residue Z is well defined (0 < Z < 1) 
by evaluating the negative slope of Im'Z,{iujn) at cj = 
[Z = (1 - ■^^^|c^=o)"^)- However, in the incoherent 
regime (T > Tcoh), the slope of Im'E{iujn) at w = 
becomes positive making the concept of Z ill defined 
(Fig. [7]). Therefore, we determined Tcoh as the temper- 
ature where the slope of the low energy self energy van- 
ishes, and found that it is almost proportional to S"^^^, 
in surprising agreement with the findings of a previous 
study of doped Mott-insulator— . 

In a two-dimensional Hubbard model, a long-range 
magnetic order at a finite T is prohibited by the Mermin- 
Wagner theorem. A FM order is possible only at T = 0. 
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At any finite T, Goldstonc modes disorder the systerrk^l^ 
and it results in a correlation length which is finite but 
exponentially large in T^^. DMFT docs not capture this 
behavior. Therefore, Tc in the context of the two dimen- 
sional model should be interpreted as an estimate of the 
temperature where the correlation length gets very large. 
In higher dimensions, wc expect a FM state at low T with 
the correct dependence on t' /t. The Nagaoka fcrromag- 
netism study using the dispersion of realistic materials 
deserves further investigations since the energy balance 
between a FM state and a PM state or the character of 
the transition is very sensitive to the details of the lattice 
structure. 



In general, n-^ — exhibits small fluctuations near the 
boundary of Tc due to the finite T. The fluctuations be- 
come especially severe through the transition from the 
FM state to the incoherent PM state near n = 1 . There- 
fore, the boundary points can be determined more pre- 
cisely by examining the temperature dependence of Xqlo 
(Fig. [5]) . Xq=o near a transition point obeys the Curie- 
Weiss form {xqlo -^T -Tc). Both coherent {n = 0.85) 
and incoherent (n = 0.95) regions show linear depen- 
dence of Xqlo ^- "^^"^ ^q=o ^ ~ ^■'^'^ barely de- 
pends on T, exhibiting Pauli paramagnetic behavior. xYol 
is greater than Xqlo ^^'^ increases as n decreases. This 
is because in DMFT xTol: ~ + Tcoh and Tcoh increases 
as n decreases^^. 



Fig. [7| shows the behavior of ImT,{iuJn) for the three 
different phases in the Tc phase diagram of Fig. [5l For 
n = 0.85 and T = 0.01, a coherent FM state is expected 
from the phase diagram. A coherent Fermi liquid is vali- 
dated by investigating the negative slope of ImT,(iLUn) at 
Lu — 0. The slope for spin a at the high frequency part is 
given by —ri-^/l — (Eq. 7) and the inequality of the 
slope indicates ^ ni confirming the FM state. The 
majority spin state has a smaller slope at high frequency 
because n_cr of the majority spin is smaller than that of 
the minority spin. Also, because the slope of the majority 
spin at w = is smaller, Z of the majority spin is larger 
than that of the minority spin. This means the quasipar- 
ticle band of the minority spin is strongly renormalized 
by correlations while the majority spin state tends to 
be similar to the non-interacting energy dispersion. For 
n = 0.85 and T ~ 0.02, a coherent PM state is estab- 
lished by observing a negative slope a,t u) — Q and no 
spin symmetry breaking. For n = 0.95 and T — 0.02, an 
incoherent PM state is expected from the positive slope 
at w = because the concept of Z is no longer valid and 
the application of Fermi liquid theory fails. Lastly, for 
fixed T = 0.02, as n increases from 0.85 to 0.95 the slope 
at high frequency also increases because increases. 




tVt 

FIG. 8; (Color online) The lowest energies of a 5=1/2 state 
and a 5=3/2 state in a f/ = oo 4-site toy model varying t' /t. 
E is the energy in units of t = 1/2. 



III. NAGAOKA FERROMAGNETISM FROM A 
4-SITE PLAQUETTE 

In order to provide a simple interpretation of why 
decreasing t' stabilizes the Nagaoka state, we examine 
the simplest possible model which retains the physics 
of the Nagaoka problem. We consider a 4-site plaque- 
tte with three electrons (one hole). The ground state of 
this model may be characterized by the quantum num- 
ber corresponding to the total spin angular momentum 
(ie. S = |, ^) and the z-direction of the spin angular 
momentum (5*2 = ±1,52 = ±^). The whole Hamilto- 
nian matrix is a 32 x 32 matrix excluding double-occupied 
sites and it is block-diagonalized to 6 distinct spin sectors 
by performing the unitary transform to the proper S, Sz 
basis. The ground state energy at each spin sector is 
determined by the exact diagonalization of Hamiltonian 
matrix. 

The lowest energy in a S" = | sector is given by —2t + t' 
and in the S = ^ sector is given by — V3t^ -I- f^. The 
energy dependence of a 5 = | state is noticeably different 
from that of a 5 = ^. In a S* = | case, doubly occupied 
states are excluded by the Pauli principle regardless of U. 
Therefore, the U = oo Hamiltonian is equivalent to the 
U — Hamiltonian where the addition of the positive 
n.n.n hopping t' contributes linearly to the increase of 
the kinetic energy. However, doubly occupied states in a 
S — ^ sector are excluded only for J7 — > oo. Therefore, 
unlike the 5 = | case, the energy dependence on t' is 
greatly reduced as the Hilbert space shrinks due to the 
infinite U. 

A 5 = I ground state is indicative of the Nagaoka fer- 
romagnetic state while a S = ^ ground state is indicative 
of a paramagnetic state. The S = ^ state is the ground 
state for t' /t < 0.24 and the energy difference increases 
approximately linearly thereafter indicating that the Na- 



7 



gaoka state is stabilized as t' /t is decreased. This is in 
qualitative agreement with the DMFT results presented 
in the previous section. The energy of the 5 = ^ state 
weakly depends on t' while the 5=1 energy decreases as 
t' /t decreases. This also explains that the stability of Na- 
gaoka ferromagnetism originates from the minimization 
of the kinetic energy. 



IV. A MEAN-FIELD SLAVE BOSON 
APPROACH 

In this section, Nagaoka ferromagnetism in a f/ = oo 
Hubbard model is studied using a mean-field slave boson 
approach. In a slave boson method, a fermion operator 
is accompanied by bosonic operators (ie. slave bosons) 
which keep track of the local occupation number. The 
three slave boson operators are e,p-\^pi and they act on 
unoccupied sites, spin-up sites, and spin-down sites, re- 
spectively. In this U = oo case, the doubly occupied sites 
are excluded. Constraints regarding the conservation of 
the occupation number are imposed with Lagrange mul- 
tipliers (A, A|, A^). The slave boson Hamiltonian is given 
by 

H = — y ^ tijcl^Zia-Z^j^Cjfj — ^ ^ ^ia{Pi(jPia ~ cJcr^'ic^) ~^ 
ij a i(7 

Y,Uplp^o + e\h~l) (10) 



where 



--e\pia 



tij=t if i,j 



are n.n, and tij=t' if i,j are n.n.n. The non- interacting 
e(k) is taken to be — 2i(cos + cos ky) — At' cos fc^ cos ky 
as in the previous section. The original Fock space has 
been enlarged including the slave boson fields. The parti- 
tion function can be calculated from the Feynman func- 
tional path integral over the original fermi fields, slave 
boson fields, and Lagrange multipliers. The integral over 
the fermi fields is straightforward because the Hamilto- 
nian is quadratic in the fermi fields. The integral over the 
slave boson fields and Lagrange multipliers should be per- 
formed using the saddle-point approximation, where the 
integral over the slave boson fields and Lagrange multi- 
pliers is approximated by putting their space and time in- 
dependent mean-field values which minimize the Hamil- 
tonian. The physical meaning of slave boson mean-field 



value is clear. The expectation value 



corresponds 



to the fraction of unoccupied sites, i.e. the hole density 
(5(1 — n). Similarly, (pjpt) equals to the spin up occu- 
pation number (rif), and {p\pi} corresponds to the spin 
down occupation number (n^). 

The free energy can be derived from the partition func- 
tion (F = —ksT In Z) and it is necessary to compare the 
free energies between ferromagnetic state and paramag- 
netic state to investigate the transition. The free energy 
is a function of magnetization m = — ni, S, and T. 



At T = 0, the free energy becomes the ground state en- 
ergy. The energies of the fully polarized ferromagnetic 
(FPFM) state (m = n^) and the paramagnetic (PM) 
state (m = 0) are given by. 



i?™/(<5) = ^6(k)e(M - e(k)) 

k 

EpMiS) = ^ E Zeik)e{p* - Zeik)) 



(11) 



(12) 



k.cr 



where Ng is the number of total sites, Z is the renor- 
malization residue given by 26/{l + 6), fi is the chem- 
ical potential in a fully polarized ferromagnetic state 
satisfying (1/A^,s) X^k ®(m ~ ^(k)) = n-f =1 — 6, and 
/i* = (/i— Ao-) is the effective chemical potential in a para- 
magnetic state satisfying {l/Ns)J2l^'^i^^* ~ ^^i^)) = 
n-f = ni = {1 - S)/2. The DOS of the FPFM state 
is the same as the non-interacting DOS (po(e)) while the 
DOS of the PM state is renormalized by a factor Z to 
1/Z-poie/Z). Unlike the DMFT method, the slave boson 
approach considers only the renormalized quasiparticle 
DOS ignoring the incoherent contribution. EpM is given 
hy Z -Eq where Eq is the non-interacting energy. In other 
words, as S reduces to 0, the energy for a paramagnetic 
state is strongly renormalized by a factor 25/{l + S) to 
avoid the doubly occupied states. That makes the FPFM 
state more stable at small 5. 

In Fig. [91(a), FPFM energy and PM energy vs n are 
shown for t' /t = 0.1, and —0.1. For all values of t' /t, the 
FPFM energy is stable at large n while the PM energy 
is stable at small n. The intermediate phase separated 
region is constructed by the Maxwell construction and 
is indicative of a fist order transition. At large n, as in 
the plaquette case, the energy curve for the paramag- 
net state depends weakly on t' while the FPFM energy 
is reduced with decreasing t' /t. This results is in qual- 
itative agreement with the previous DMFT results. As 
t' /t decreases, the FPFM state becomes more stable and 
the critical density, Uc decreases. Just as in the DMFT, 
the large spectral weight of the non-interacting DOS at a 
low energy makes FPFM state energetically favorable at 
t' /t = —0.1. When t' is 0, the energy difference between 
FPFM and PM vanishes at ric=2/3, in agreement with 
the previous slave boson calculations^. 

We also calculate the inverse of uniform magnetic sus- 
ceptibility (x~^) to study the instability of the PM state. 
The analytic expression is 



1 



X \m,=0 



-y- 

Ns V (1 



[i + sy 



2p(/i*) 1 + 5 
rZe(k)e(/i*-Ze(k)) 



(13) 



where p{^J■*) is the renormalized DOS given by 1/Z ■ 

The trends in x -'^ are consistent with the results shown 
in Fig. [ni(a). As t' /t decreases, spin susceptibility di- 
verges at smaller density (see Fig. [HI (b)). However, the 
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n 

FIG. 9: (Color online) (a) Fully polarized ferromagnetic 
(FPFM) energy and paramagnetic (PM) energy vs n varying 
il jt (0.1 (top), (middle), and -0.1 (bottom)) Inset : Maxwell 
construction to determine the PS region, (b) The inverse of 
the uniform magnetic susceptibility (x^^) a-t 'm=(i vs n vary- 
ing t' jt (0.1, 0, and -0.1). (c) The chemical potential (^) vs 
n at i! jt = 0.1, 0, and -0.1. 



divergence of the spin susceptibility does not coincide 
with the thermodynamic phase transition identified by 
the total energy differences. The phase transition is thus 
always first order within the slave boson approach. 







i! lt=-^.\ 


t'/t^O 


t7t=0.1 


DMFT 
(T = 0.01) 


Tic 


0.705 


0.815 


N/A 


order 


First 


Second 


N/A 


Slave boson 

(r = o) 




0.53 


0.67 


0.83 


order 


First 


First 


First 



TABLE L ric and the order of the ferromagnetism transition 
in a (7 = oo Hubbard model from both the DMFT-(-CTQMC 
approach and the slave boson approach with t' /t= -0.1, 0, 
and 0.1. N/A means no transition to FM state occurs. 



Fig. [H] (c) shows that a flat chemical potential region 
exists at any t' /ts in a /u vs n graph. This is a generic fea- 
ture of a first order transition and this region represents 
the coexistence of the FPFM and PM phase. This coexis- 
tence region is larger for negative t' /t favoring transition 
to the FPFM phase. 



V. COMPARISON OF THE SLAVE BOSON 
RESULT AND THE DMFT+CTQMC RESULT 




0.6 0.7 0.8 0.9 1 



n 

FIG. 10: (Color online) Paramagnetic energy from both the 
DMFT-HCTQMC (T = 0.01) and the slave boson (T = 0) 
approach vs n at t' /t=-Q.l (the top panel), (the middle 
panel), and 0.1 (the bottom panel). 
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0.65 0.7 ^ 0.75 0.8 

FIG. 11: (Color online) The renormalization residue (Z) of 
the slave boson method and the DMFT+CTQMC method 
it' = 0). 



The slave boson method overestimates the region of 
the stable FM state as compared to DMFT and it fa- 
vors a first order transition (see Table 1). This is be- 
cause the slave boson approach overestimates the para- 
magnetic kinetic energy as compared to the DMFT ap- 
proach (Fig. [TUl) . The quasiparticle residue Z of the 
DMFT approach is evaluated by (1 - -^^t^U^o)^^ on 
the imaginary frequency axis while Z of the slave bo- 
son approach is given by 26/(1 + 5). Fig. [Til shows that 
Z of the slave boson study is overestimated as compared 
to the DMFT+CTQMC case. The slave boson technique 
used in this paper is based on the mean-field saddle-point 
approximation and it does not treat the strong corre- 
lation effect properly. Even though DMFT ignores the 
spatial correlation effect, the temporal correlations arc 
treated exactly by CTQMC. Moreover, the mean-field 
slave boson approach evaluates the total energy as the 
sum of coherent quasiparticle energies (Eq. [12]) while the 
total energy of DMFT-hCTQMC includes contributions 
from both the incoherent and coherent effects. The over- 
estimated Z in the slave boson case underestimates the 
kinetic energy while the ignorance of contribution from 
the incoherent part overestimates the energy. As a result, 
the two errors of the slave boson approach cancel each 
other giving a slightly overestimated energy as compared 
to the DMFT+CTQMC result. 

Additionally, the graph in the slave boson method 
almost coincides with the DMFT+CTQMC result com- 
paring Fig. [H and Fig. [5| (b). It is not certain 
how the renormalization residue Z affects iii the 
DMFT+CTQMC case, and the contribution from the 
incoherent part is also unclear. Therefore, further study 
will be required to fully understand the positive agree- 
ment of X ill the two methods. 



VI. CONCLUSION 

To summarize, we investigated Nagaoka ferromag- 
netism in the U = oo Hubbard model including n.n hop- 
ping t and n.n.n hopping f . This model was solved using 
DMFT with CTQMC, and the mean-field slave boson 
approach. Even a small value of t' /t yields a significant 
impact on the stability of Nagaoka ferromagnetism. The 
DMFT results show that the FM state is more stable for 
negative t' /t, and this is supported by the slave boson 
method (see Tabic 1) and can also be understood from 
diagonalization of the 4-site plaquette. The energy of the 
minimum spin state {S = 1/2) depends weakly on t' /t, 
while the energy of the maximum spin state (5" = 3/2) 
depends linearly on t' /t. Therefore, the maximum spin 
state becomes more stable for negative t' /t. 

In both slave boson and DMFT methods, the para- 
magnetic energy does not vary much with t' /t due to 
the strong renormalization of the quasiparticle band 
(see Fig. [TU|) . However, the fully polarized ferromag- 
netic energy depends on t' /t in a similar fashion as the 
non-interacting kinetic energy since the correlations are 
weaker in the broken symmetry state. The negative t' /t 
gives a high spectral peak in the occupied part of the 
spectra of the non-interacting system. As a result, the 
energy of the FM state is lower and the ferromagnetism 
is stabilized in this case. 

Within DMFT, the nature of the transition also varies 
with t' /t. A first order transition accompanied by the PS 
of the FM and PM state occurs at t' /t = —0.1 while a 
second order transition occurs at t' /t ~ 0. In the slave 
boson approach, the transition is always first order re- 
gardless of t' /t. This is because the slave boson method 
overestimates the PM energy. The DMFT result shows 
that when n — > 1, the FM state becomes unstable as T 
exceeds Tcoh- In other words, ferromagnetic state is only 
stable within the coherent Fermi liquid regime. 

The U = oo one band Hubbard model is a toy model 
and does not describe any specific material. However it 
is physically realizable in an optical lattice, due to the 
recent developments in controlling cold atoms in optical 
trap o^^i^^ . These systems are highly tunable, and the 
hopping parameter t and the on-site interaction U can 
be adjusted by varying the ratio of the potential depth of 
the optical lattice to the recoil energy (Vo/Eji) or the ra- 
tio of interatomic scattering length to the lattice spacing 
(as/d). In order to realize the one-band Hubbard model 
with a large U {U/t > 100), Vq/Er « 30 and Us/d < 0.01 
should be the range of parameters in the optical lattice 
(See Fig. 4 of Ref. 29). The tuning of the next-nearest 
neighbor hopping t' can be achieved by engineering opti- 
cal lattices with a non-separable laser potential over each 
coordinate axis. 

It will be very interesting to test these DMFT results 
experimentally. Usually, the atomic trap potential is ap- 
plied to confine atoms in the optical lattice, and the po- 
tential varies smoothly having the minimum at the center 
of the trap. The phase separation between the FM and 
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the PM phase at t' /t = —0.1 (taking place between the 
densities n = 0.696 — 0.715) can be observed in the op- 
tical lattice as three spatially separated distinct regions. 
The atom-rich FM region will tend to move to the center 
of the optical lattice to be energetically stabilized while 
the hole-rich PM region will reside on the edge of the 
optical lattice. Since the total spin is a conserved quan- 
tity, the FM region will be located at the center of the 
trap and will consist of two domains containing the up 
or down species. Raising the temperature will destroy 
the ferromagnetic magnetic state and consequently the 



spatial patterns within the trap. 
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